/*!  \file  BiSVD.h
 *	  \brief  
 * 
 *   \details   details
 *   \author    Yordan Madzhunkov German Research School for Simulation Sciences GmbH
 *   \version   1.0
 *   \date      07/02/2013 04:23:03 PM
 *   \copyright GNU Public License.
 */
#pragma once

#include "ImplementationFactory.h"
#include "UnitTester.h"
class TriEigen;
class BiSVD : public Implementation { 					// interface for bidiagonal SVD
public:
	// Computes the SVD of bi diagonal matrix B 
	// B = U2 * Sigma * V1
	// B is a Bidiagonal Matrix
	// U2, V2^T are orthogonal matrices
	void saveCompute(const Matrix & B, Matrix & U, Matrix & Sigma, Matrix & V); 
	virtual ~BiSVD() {};
	BiSVD(const char * name_);
	virtual bool setNumberOfThreads(int nThreads_);
protected:
	void setTriEigen(TriEigen & implementation);
	virtual void compute(const Matrix & B, Matrix & U2, Matrix & Sigma, Matrix & V2T ) = 0;
	void computeNormal (const Matrix & B, Matrix & U, Matrix & S, Matrix & VT);
	void computeNormal2(const Matrix & B, Matrix & U, Matrix & S, Matrix & VT);
	void computeGK     (const Matrix & B, Matrix & U, Matrix & S, Matrix & VT);
	bool checkInput(const Matrix & A, Matrix & U,
 	 					 const Matrix & S, Matrix & VT);
	TriEigen * eigen;
	bool aUpper, aLower, aURowMajor, aUColMajor, aVRowMajor, aVColMajor;
	bool VColVec; // does computation stores vectors in rows or in cols
};
#define BISVD(x) \
class x : public BiSVD { \
public: \
	x();\
	virtual ~x(){}; \
	virtual void compute(const Matrix & B, Matrix & U2, Matrix & Sigma, Matrix & V2T ) ; \
};
#define BISVD_NORMAL(x) \
class x : public BiSVD { \
public: \
	x();\
	virtual ~x(){ if (eigen) delete eigen; eigen = 0;}; \
	virtual void compute(const Matrix & B, Matrix & U2, Matrix & Sigma, Matrix & V2T ) { computeNormal(B, U2, Sigma, V2T);} \
};

#define BISVD_NORMAL2(x) \
class x : public BiSVD { \
public: \
	x();\
	virtual ~x(){ if (eigen) delete eigen; eigen = 0;}; \
	virtual void compute(const Matrix & B, Matrix & U2, Matrix & Sigma, Matrix & V2T ) { computeNormal2(B, U2, Sigma, V2T);} \
};
#define BISVD_GK(x) \
class x : public BiSVD { \
public: \
	x();\
	virtual ~x(){ if (eigen) delete eigen; eigen = 0;}; \
	virtual void compute(const Matrix & B, Matrix & U2, Matrix & Sigma, Matrix & V2T ) { computeGK(B, U2, Sigma, V2T);} \
};

// BENCHMARKS
template <typename T> class Metric;
class MatrixGenerator;
class BiSVDBenchmark : public Benchmark {
public:
	BiSVDBenchmark(BiSVD & bisvd_);
	virtual ~BiSVDBenchmark();
	virtual void setRunConfiguration(const int n) ;
	virtual int  getProblemSize(const int configuration) const;
	virtual void prepareRun();
	virtual void makeRun();
	virtual void finalizeRun();
	enum { nSvdNorms = 2 };
protected:
	void addSVDNormsMetrics();
	void measureSVDNorms();
	DISALLOW_COPY_AND_ASSIGN(BiSVDBenchmark);
	MatrixGenerator  * generator;
	BiSVD            * bisvd;
	Matrix * Bi, * U, * V, *S;
	Metric<double >  * svdNormMetrics[2];
};

class BiSVDScaling : public Benchmark {
public:
	BiSVDScaling(BiSVD & bisvd_);
	virtual ~BiSVDScaling();
	virtual void setRunConfiguration(const int n);
	virtual int  getProblemSize(const int configuration) const;
	virtual void prepareRun();
	virtual void makeRun();
	virtual void finalizeRun();
	enum { nSvdNorms = 2 };
protected:
	void addSVDNormsMetrics();
	void measureSVDNorms();
	DISALLOW_COPY_AND_ASSIGN(BiSVDScaling);
	MatrixGenerator  * generator;
	BiSVD            * bisvd;
	Matrix           * Bi, * B2, * U, * V, *S;
	Metric<double >  * svdNormMetrics[2];
};

// UNIT TESTS
class BiSVDUnitTest : public UnitTest { 
public: 
	BiSVDUnitTest ( BiSVD &); 
	virtual ~BiSVDUnitTest (){}; 
	virtual bool test() const; 
private:
	bool test1(const int m, const int n, const double * d, const double limit,
		const double * s_ref = NULL) const;
	DISALLOW_COPY_AND_ASSIGN(BiSVDUnitTest);
	BiSVD & bisvd;
};
